from pyproj import Transformer, CRS

# 定义CGCS2000和WGS84的坐标系
crs_cgcs2000 = CRS.from_epsg(4490)  # CGCS2000坐标系
crs_wgs84 = CRS.from_epsg(4326)    # WGS84坐标系

# 创建转换器
transformer = Transformer.from_crs(crs_cgcs2000, crs_wgs84)

# 示例CGCS2000坐标点 (x, y, z)
cgcs2000_points = [
    (773118.408,3998395.095,292.133),
    (773118.393,3998394.896,292.344),
    # 添加更多点...
]


# 转换坐标点
wgs84_points = []
for x, y, z in cgcs2000_points:
    lat, lon, alt = transformer.transform(x, y, z)
    wgs84_points.append((lat, lon, alt))

# 输出转换后的WGS84坐标点
for i, (lat, lon, alt) in enumerate(wgs84_points):
    print(f"Point {i+1}: Latitude={lat}, Longitude={lon}, Altitude={alt}")

# 如果需要将WGS84坐标转换为平面坐标系（如UTM），可以继续转换
# 例如，转换为UTM Zone 50N
crs_utm = CRS.from_epsg(32650)  # UTM Zone 50N
transformer_utm = Transformer.from_crs(crs_wgs84, crs_utm)

utm_points = []
for lat, lon, alt in wgs84_points:
    easting, northing = transformer_utm.transform(lat, lon)
    utm_points.append((easting, northing, alt))

# 输出转换后的UTM坐标点
for i, (easting, northing, alt) in enumerate(utm_points):
    print(f"Point {i+1}: Easting={easting}, Northing={northing}, Altitude={alt}")